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ABSTRACT 

We present a new formalism to describe the ratios and profiles of emission 
lines from hydrogen in Balmer-dominated shocks. We use this model to interpret 
the measured widths and ratios of broad and narrow Ha, B.f5 and Lya emission 
lines in supernova remnants (SNRs). Our model results agree fairly well with 
those obtained previously by Chevalier, Kirshner & Raymond (1980) and are 
consistent with observations of several SNRs. The same model fails to account 
for the ratio of broad to narrow line emission from the reverse shock in SNR 
1987A as observed by Heng et al. (2006). We suggest that this discrepancy 
between theory and observation results from a faulty assumption that Balmer- 
dominated shocks can be treated as sharp discontinuities. If the spatial structure 
of the shock transition zone is taken into account, the predicted ratios of broad to 
narrow line emission in most SNRs will change by modest factors, but the ratio 
in SNR 1987A will increase substantially. Significantly greater shock velocities 
will be required to account for the observed full widths at half-maximum of the 
broad emission lines in most SNRs. 

Subject headings: methods: shock waves — supernova remnants 



1. INTRODUCTION 



When a fast astrophysical shock enters neutral interstellar gas, we see an optical spec- 
trum that is dominated by Ha and other Balmer lines and is nearly devoid of the forbidden 
lines typically seen in emission line nebulae. The Balmer emission is caused by impact exci- 
tation of neutral hydrogen atoms by fast ions and electrons in the shocked gas. The emission 
lines in such "Balmer-dominated" or "non-radiative" shocks have two components. The first 
is a narrow component with a line width characteristic of the cold interstellar gas, which re- 
sults from direct excitation of the neutral hydrogen atoms. The second is a broad component 
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with a line width characteristic of the thermal velocity broadening of the shocked protons. 
The broad component results from charge transfer reactions, which pass electrons from the 
nearly stationary hydrogen atoms to the shocked protons. The resulting fast hydrogen atoms 
may be created in excited states or may be excited by subsequent collisions with other fast 
ions and electrons. The resulting ratio of the broad to narrow components and the profiles 
of the broad Balmer lines depend on the equilibration of electron and ion temperatures in 
the shocked gas, on resonance trapping of the Lyman lines, and on the velocity and the 
inclination of the shock surface. 

The theory to interpret the emission spectra of such Balmcr-dominatcd shocks was de- 
veloped originally by Chevalier & Raymond (1978, hereafter CR78) and Chevalier, Kirshner 
& Raymond (1980, hereafter CKR80) and augmented by Smith et al. (1991, hereafter S91) 
and Ghavamian et al. (2001, hereafter GOl). It has been used to interpret the emission 
spectra of several supernova remnants (SNRs), including: SN 1006 (CKR80; Kirshner, Win- 
kler & Chevalier 1987, hereafter KWC87; S91; Winkler & Long 1997; Ghavamian et al. 
2002, hereafter G02); Kepler (Fesen et al. 1989, hereafter F89; Blair, Long & Vancura 1991, 
hereafter BLV91; Sollerman et al. 2003; Sankrit et al. 2005); Tycho/SN1572 (KWC87; S91; 
Ghavamian et al. 2000; GOl), RCW 86 (Long & Blair 1990, Smith 1997; GOl; SOS), portions 
of the Cygnus Loop (Raymond et al. 1983; Fesen & Itoh 1985; Hester, Raymond & Blair 
1994; GOl), and four remnants in the Large Magellanic Cloud (Tuohy et al. 1982; S91; Smith 
et al. 1994). Generally, the theoretical models fit the observations fairly well, but significant 
discrepancies (by factors ~ 2) between the theoretical and observed broad-to-narrow line 
ratios persist in some cases, as we shall discuss. 

Supernova 1987A has provided a dramatic new example of a Balmer-dominated shock. 
In contrast to the SNRs cited above, broad Balmer emission from SNR 1987A comes from 
hydrogen atoms in the rapidly expanding (~ 12,000 km s~^) supernova debris rather than 
the stationary interstellar gas. In this case, the highest velocity component of the emission 
lines, called "surface emission," results from the excitation of freely-streaming hydrogen 
atoms that cross the reverse shock. Observations of the profiles of Ha and Lya emission 
from the reverse shock with the Space Telescope Imaging Spectrograph (STIS) have enabled 
us to measure the rapidly evolving fiux of hydrogen atoms across the reverse shock and to 
reconstruct partially the three-dimensional geometry of the shock surface (Michael et al. 
1998a, b, 2003, with the last hereafter M03; Heng et al. 2006, hereafter H06). In addition to 
the surface emission, H06 identified a lower velocity component of the Ho; and Lya emission 
lines, which they called "interior emission." Evidently, this interior emission results from 
charge transfer of electrons from the freely-streaming hydrogen atoms to protons in the 
shocked gas and is analogous to the broad component seen in the Balmer lines of other 
SNRs. We also note that Balmer-dominated bow-shocks have been observed in four pulsar 
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nebulae (Bucciantini 2002). 

Here, we revisit the theory of B aimer-dominated shocks. We hope to resolve the dis- 
crepancies between observations and theory, to understand better how to interpret the line 
profiles, and to extend the existing theoretical models to the very high shock velocities seen 
in SNR 1987A. In §2, we state the assumptions of our model and mention their limitations. 

In §3, we discuss the cross sections used and present the resulting velocity distribution func- 
tions and rate coefficients. In §4, we present a new formalism for deriving the intensities and 
the profiles of broad and narrow emission lines. In §5, we compare our results to previous 
models and to existing data on Balmer-dominated remnants and SNR 1987A. Finally, in §6, 
we discuss the limitations of the model and future development that is needed. 



2. MODEL AND ASSUMPTIONS 

The basic model we consider is a plane-parallel shock with velocity Vg that strikes a 
stationary gas of cold neutral hydrogen and helium, with the latter having fractional abun- 
dance XHe = nHe/nn- (Other elements that may be present with typical cosmic abundances 
will have no significant effects.) We assume that the shocked gas is fully ionized and that 
the ions and electrons comprise a fluid having a velocity relative to the shock frame. 
We assume that the ions and electrons have Maxwellian velocity distributions. 

There is an uncertainty regarding the equilibration of post-shock ion temperatures. 
Electrons and protons are heated to temperatures having a minimum ratio equal to that of 
their masses (i.e., Tg/Tp ~ 2000). The Coulomb equihbration time often exceeds the age 
of the remnant, in which case the temperatures will remain unequal throughout the transi- 
tion zone where the neutral atoms become ionized. However, plasma waves and magneto- 
hydrodynamical turbulence at the shock front may transfer energy rapidly from protons to 
electrons. We consider two models. In the first, "Model F", the electron and ion temper- 
atures are fully equilibrated. In the second, "Model N", two-stream plasma instabilities 
produced by protons refiected upstream result in ~ 0.25Tp (Cargill & Papapoulos 1988; 
see also references in G02). Following G02, we define an equilibration parameter, feq, such 
that 
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where ji^ ^XHe) /(2 + 3xije) is the mean molecular weight. For example, ji 0.61 for 
Galactic abundances {xhc — 0.1). Consequently, f^q — 1 for Model F and feq ~ 0.35 for 
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Model N. We further assume that equihbration between the protons and the alpha particles 
occurs so rapidly that their post-shock temperatures are equal, namely Tp — T^- 

The most questionable assumption of this model is that the hydrogen atoms enter a 
shocked ion plasma having uniform density, velocity, and temperature given by the adiabatic 
jump conditions. In reality, the hydrogen atoms are ionized and excited in a transition zone 
where the ions are produced, decelerated, and heated. No model, including the present 
one, of emission from Balmer-dominated shocks has included a proper description of this 
transition zone. In §6, wc discuss the consequences of this approximation further, deferring 
a more detailed treatment to a future paper. 



3. REACTION RATES AND CROSS-SECTIONS 

3.1. FREELY-STREAMING ATOMS 

Consider a beam of freely-streaming, neutral hydrogen atoms in the ground state that 
crosses the shock and encounters a thermal plasma (with relative velocity Vgh = 3fs/4). Such 
an atom has three possibilities open to it — electron or ion impact excitation (denoted ) , 
impact ionization (denoted "/"), and charge transfer (denoted "T"'). Unless specified, when 
the label "T" is used it refers to charge transfer to both the ground and the excited states, 
while "T*" denotes only the latter. 

The rate that a hydrogen atom will have a reaction X (where X stands for E, 7, or T) 
with a particle of type s (in units of s""-*^) is given by 



R{Vsh; CTx) =na fa{Va) fb{Vb) CTxiWa - Vb\) \Va - Vb\ d^Va Svj,. (2) 




In the case at hand, the ions (with the subscript s denoting protons, electrons or al- 
pha particles) are assumed to have mass m^, number density n^, temperature T,, and a 
Maxwellian velocity distribution function, 

(Tfl \ 

where /m,o = {'ms/27rkTs)^^'^. The beam of hydrogen atoms has the distribution function, 

fo{vH, Vsh) = S{vh - Vsh)- (4) 

The rate of the reaction X between them is then 

Rxoivsh] CTx) = Usfufi J J exp (^--^^^ ^i^H - Vsh) crx{\vH - v\) \vh - v\ d^VH d^v 
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in cylindrical coordinates, where d^v — 2'KVr dvr dvz and Av = v"^ + [v^ — VghY- The 
reaction rate coefficient is simply Rxq = Rxo/ng- We further define the fractional abundance, 
Xs = ris/riH, such that Rxo/riH = XsRxq- 

A fraction of the incoming beam of hydrogen atoms will undergo charge transfer reac- 
tions with protons in the shocked plasma, producing hydrogen atoms having a distribution 
function given by 

fl{V:Vsh) = Si{Vsh) fmiv) \V-Vsh\ '7T{\V-Vsh\), (6) 

where the charge transfer cross section, (Jr, includes reactions to all excited states as well 
as to the ground state. The value of the normalization constant, Si{vsh), is fixed by the 
condition J fi d^v — 1. If we define 

RToiVsh^CTT) = J fM{v) \V-Vsh\ Ot{\v - Vsh\) dh, (7) 

then it follows that Si — R^^. 

The function fi{v,Vsh) describes the velocity distribution function of hydrogen atoms 
resulting from charge transfer reactions of the original hydrogen beam with protons in the 
ionized gas. The subscript "1" denotes that the atoms are the result of one charge transfer 
reaction. The distribution fi is not Maxwcllian because it includes factors of the relative 
velocity between the beam and the ions and the charge transfer cross section. For example, 
it has the value fi{vsh, Vsh) = at the velocity of the initial beam. 



2kT, 



Av ax{Av) v,. dv,. dv^, (5) 



3.2. SUBSEQUENT REACTIONS 

These new hydrogen atoms can in turn be excited, ionized, or undergo another charge 
transfer reaction. If they are excited, they will emit one or more photons and remain in the 
same distribution function. 

The velocity distribution function of hydrogen atoms that have undergone i charge 
transfer reactions is given by 

/.(«>-..) = /„(«)/ .AW- ff\) dV. (8) 

where — RtI_-^ ■ The corresponding rate coefficients for reactions X of such atoms are 

Rxiivsh,(^x) =11 fuiv) fi{vH,Vsh) \vh -v\ ax{\vH -v\) d^VH d^v. (9) 
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3.3. 



CROSS SECTIONS 



To calculate the reaction rate coefficients for excitation, ionization, and charge trans- 
fer (by interactions with protons, electrons, and alpha particles), we use analytical fitting 
functions to the cross sections provided by Janev & Smith (1993). 

We must distinguish between the net charge transfer cross section, ctt, and the cross 
sections, (Jt*, for charge transfer reactions which leave the hydrogen atom in excited states 
(and consequently produce photons). For the latter, we use the cross sections of Barnett 
(1990), who tabulates results for n — 2s and n — 2p. To obtain cross sections for n > 2, we 
use the scafing relation (Janev & Smith 1993; CD. Lin 2006, private communication). 



We include all excitation and charge transfer reactions up to and including n — 6. 

We also consider charge transfer reactions between alpha particles and hydrogen atoms, 
which produce singly-charged helium ions. We include reactions to He"*" in the states n — 
Is, 2s, 2p, 3s, 3p, 3d and 4p, using the data of Barnett (1990). We do not follow the evolution 
of these ions; instead, we regard the reactions producing them as equivalent to impact 
ionization of hydrogen and include the corresponding cross sections as a contribution to the 
net impact ionization cross section. 

Examples of these cross sections are shown in Figure 1. 



Figure 2 shows examples of one-dimensional velocity distribution functions of hydrogen 
atoms resulting from charge transfer reactions for two representative shock velocities. The 
velocity axes represent the component of velocity normal to the shock surface; the distribu- 
tion functions have been integrated over the transverse velocity components. For each case, 
the symmetric, thin, solid curve represents the Maxwellian distribution of protons in the 
shocked gas, while the dotted and/or dashed curves represent the one-dimensional hydrogen 
atom velocity distribution functions, ^^(fx) = / fi{v) dvxdvy, that result from the ith charge 
transfer reactions. The thick, solid curve represents the projected composite distribution of 
all hydrogen atoms in the shocked gas, 4'comp-i resulting from one or more charge transfer 
reactions (to be discussed in §4). The two lower panels represent Model F, in which the 
electron and ion temperatures are fully equilibrated, while the upper panels represent Model 
N, in which = Q.2bTp (§2). 
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3.4. 



DISTRIBUTION FUNCTIONS AND RATE COEFFICIENTS 



For shock velocities Vg < 2000 km s~^, the atomic distribution functions, /j, resulting 
from charge transfer reactions differ only slightly from the Maxwellian distributions, /m, of 
the protons. We call these /j "skewed" Maxwellian distributions. However, for velocities 
Vg > 2000 km s~^, the /i have peaks that are shifted in velocity toward the original beam, 
progressively more so as the shock velocity increases. This behavior is a consequence of the 
fact that the charge transfer cross section decreases rapidly for relative velocities > 2000 km 
s~^ (Fig. 1). With each subsequent charge transfer, the /j shift towards /m- 

Given the distribution functions and cross sections, we calculate the rate coefficients 
from equations (5) and (9). These coefficients arc displayed as functions of the shock velocity 
in Fig. 3, for both the original beam and the atoms that result from the ffist charge transfer 
reaction. Rate coefficients involving atoms that have undergone two or more charge transfers 
are nearly indistinguishable from those resulting from the first charge transfer. Note that for 
shock velocities 200 ^ Vg 2000 km s~^, the rate coefficients for excitation, ionization and 
charge transfer are weakly dependent on shock velocity and are comparable in magnitude. 
Note also that the excitation and ionization rate coefficients are dominated by electron 
collisions. However, for shock velocities Vg ^ 2000 km s~^ the net rate coefficients for 
excitation and ionization continue to increase, while the rates for charge transfer decrease 
rapidly. 



4. THE REACTION TREE 

4.1. REACTIONS PER ATOM 

To calculate the ratios and profiles of hydrogen lines from a Balmer-dominated shock, 
we must account for all possible interactions of the hydrogen atoms. We wish to track the 
fate of every hydrogen atom that crosses the shock, and to calculate the number and profile 
of every photon produced by the atom. The possibilities are illustrated by the "reaction 
tree" shown in Fig. 4. 

Consider a hydrogen atom belonging to a distribution function, /j, resulting from i 
charge transfer reactions (or "skewings" ) . The rate at which this atom will have an excitation 
(Ei) or ionization (/j) is given by the sum of the rate coefficients over species (protons, 
electrons and alpha particles). The rate at which the atom will have a charge transfer (Tj) 
reaction is given by Rt^ = npRT^ ■ Therefore, the probability that this atom will have another 
charge transfer reaction before it experiences an excitation or ionization is given by 
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where i?^;. is the net rate of excitations and i?/. is the net rate of ionizations. 

The net rate at which an atom will have an impact excitation after the ith charge 
transfer reaction is given by 

REi — ^ |^^e-REi,p,n + npREi,e,n + naREi,a,n , (12) 
n 

where REi,n is the rate coefficient for excitation to the level n by protons {p), electrons (e) 
or alpha particles {a) (Fig. 3). The net probability that the atom will be excited to any 
level n is given by 

Rsi + Rii + Rt, ^ ^ 

Likewise, we define the probability, P/., that the atom will have a coUisional ionization after 
the ith charge transfer reaction. 

Figure 5 displays these probabilities for i = and i = 1 as functions of the shock 
velocity. The probabihties for i > 1 are almost identical to those for i = 1. In practice, 
we average over the rate coefficients from i = 1 to i = 4 to obtain mean rate coefficients, 
from which wc construct the probabilities. Note that for shock velocities Vg ^ 1000 km 
s~^, > Pe, > Pii- These probabilities imply that every hydrogen atom that crosses a 
shock with Vg ^ 1000 km s~^ will most likely have one or more charge transfer reactions and 
excitations before it is ionized. 

With these probabilities in hand, wc can track a distribution of hydrogen atoms as it 
proceeds through the reaction tree (Fig. 4). First, consider the average number, Nt^., of 
charge transfers that a hydrogen atom in the original beam will have. Before doing so, it 
can have any number of excitations, so 

oo 

Nt, = Pto E(^^o)' = nrV- (14) 

i=0 ^0 

The average number of ionizations by the atom in the original beam before it has one charge 
transfer is given by 

Pto Pio 



Likewise, for excitations. 



J- - -TEo ~ ^Eo 



Neo = E(^^o)^ = (16) 



After the first charge transfer, the probabilities, Pxf, for excitation, ionization, and 
charge transfer are nearly independent of the number of charge transfers, so we may drop 
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the subscript i to describe the number of reactions per atom for the remainder of the tree. 
(However, see Appendix A.l for the exact expressions.) So, the average number of atoms 
which will survive to the ith charge transfer reaction is given by 



To 



l-Pf 



1 - -Peo 

The average number of charge transfer reactions that will occur is given by 

-1 



To 



l-Pf 



E 



^0 i=l 

The average number of excitations is 



E 



1-P, 



Eo 



Peo + 



Pi 



PtqPe 
Pi 



E 



Eo 



(17) 
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(19) 



Naturally, Nj — 1. The number of reactions per atom, for I and T, are plotted as a 
function of the shock velocity in Fig. 6. 



4.2. LINE STRENGTHS AND PROFILES 

To find the intensities and profiles of emission lines, we multiply the rate at which atoms 
are excited by the appropriate weighting factor that a given line photon will be emitted after 
that reaction. For example, for narrow line emission we find a rate per atom of 

4(Ha) = — ^ V REo,n C^s (20) 

and 

^ m 

4(Lya) = — J2 REo,n Cn2, (21) 

where the C„„/ are the cascade matrices (Seaton 1959) representing the probabihty that 
a hydrogen atom excited to state n will make a transition to the state n' < n via all 
cascade routes. In Appendix A. 2, we present analytical expressions for C„„/ and note that 
excitations are not necessarily distributed statistically among the angular momentum states. 
The quantity Req^u is the rate coefficient for excitations to the state n before any charge 
transfer reaction occurs. 

The broad line emission has contributions from all excitations after the first charge 
transfer and also from charge transfer reactions to excited states: 



„ m 



{RE,n + RT*,n) + Rt*, 



Cn3 (22) 
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and 



1 



m 




l-P, 



E 



{RE,n + P'T*,n) + ,n C, 



(23) 



n=2 



where RE,n is the rate of excitations to state n after any number i > 1 of charge transfer 
reactions and i?T*,n is the corresponding rate of charge transfer reactions to excited states. 

We construct a composite distribution function, fcomp, for the hydrogen atom velocity 
distribution by summing the fi, weighted by the probabihties that the atoms undergo i 
charge transfers: 



Figure 2 shows the function (t)comp, which is fcomp integrated over velocities transverse to the 
shock, for two representative shock velocities. Note that for Vg = 500 km s^^, (f)cornp differs 
only slightly from the projected proton distribution, (pM, hence validating the assumption 
made in earlier models. In contrast, for Vg — 10,000 km s~^, (pcomp differs greatly from 
0M- Its peak is shifted substantially toward the velocity of the original beam. In fact, it is 
dominated by 0i. It contains a small contribution from 02 and negligible contributions from 
<f>i>2, owing to the small probabilities that the atoms will undergo multiple charge transfer 
reactions. 

For a face-on shock, the velocity profile of the broad line is given by (pcomp- Since fcomp is 
not isotropic, the line profile will depend on the aspect (or viewing) angle, 9y. To investigate 
this angle dependence, we perform a coordinate transformation from cylindrical to Cartesian 
coordinates and then rotate the distribution about the x-axis by 9^- Figure 7 shows the full 
width at half-maximum (FWHM) of the broad emission line as a function of shock velocity 
for shocks viewed face-on {6^ = 0) and edge-wise {6^ = 90°). Because the ions are hotter 
in the non-equilibrated Model N than in the fully-equilibrated Model F, we see that the 
FWHM is greater (by a factor of between 1.1 and 1.2) in model N than in Model F. We also 
see that the FWHM increases almost hnearly with shock velocity for Vg < 2000 km s~^ but 
that it begins to level off for Vg > 2000 km s~^, consistent with earlier studies. In general, 
the FWHM is a weak function of 9^; it is shghtly greater (by a factor 1.1) for edge- wise 
than for face-on shocks. 




(24) 



5. 



RESULTS 



The observable quantities for Balmer-dominated shocks are the profiles (or widths) of 
the broad lines, the velocity shifts between the peaks of the broad line and narrow lines, and 
the intensities of the broad and narrow lines. Below, in §5.1, we discuss the interpretation 
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of these quantities for several well-known SNRs in the light of the models described above 
and we compare our results to previous work. Then, in §5.2, we discuss how our results may 
be used to interpret the emission from the reverse shock in SN 1987A. 

5.1. BALMER-DOMINATED REMNANTS 

As CKR80 and others have emphasized, we can infer the velocity, Vg, of a Balmer- 
dominated shock from the observed FWHM. Our results are fairly close to previous ones. 
For example, for our fully equilibrated Model F (viewed edge- wise), and Vg — 2000 km s~^, 
we find FWHM = 1600 km s~^ (Fig- 7), which may be compared to corresponding values 
in the literature: FWHM ^ 1800 km s"^ estimated from Fig. 3 of CKR80, FWHM ^ 1600 
km s^^ estimated from Fig. 7 of S91, and FWHM ^ 1400 km s^^ estimated from Fig. 9 of 
GOl. Unfortunately, we have not been able identify the sources for the differences among 
these values. However, we note that there is a typographical error in equation (5) of CKR80, 
which propagated into papers hke KWC87; this was later corrected in S91. 

In Table 1 (column 3), we list values of the FWHM of the broad Ha line measured in 
several SNRs. Note that values measured by different authors for a given SNR do not neces- 
sarily apply to the same part of the shock front, so the differences among such measurements 
may be real. We also list in column 4 the range of shock velocities inferred by those authors 
by fitting their models to the observations. This range takes into account uncertainties both 
in the measured FWHM and in the model itself, with the lower values of Vg coming from 
non-equilibrated models and the higher values from fully equilibrated ones. In columns 6 
and 7, we list the values of Vg that we infer by fitting our models to the same data. The 
values in column 6 are inferred from the non-equilibrated Model N, while those in column 
7 are from the fully equilibrated Model F, and the errors in each case derive only from the 
errors in the measured values of the FWHM. (In making such fits, we used the model curve 
for edge-wise aspect, which we beheve to be close to the actual situation in most cases.) In 
general, we see that the values of Vg inferred from our models are in fairly good agreement 
with those inferred by other authors, but that in some cases we infer values of Vg that are 
lower than the original values by as much as 28%. 

The velocity offset, Av{vs, 6^), between the centroids of the narrow and broad compo- 
nents is sensitive both to the shock velocity, Vg, and the aspect angle, 9^, of the shock. We 
can write 

A^;(^;„ 9^) = A^;(^;„ 0) cos(^^), (25) 

where Ai>(vs,0), the offset for a face-on shock, is shown in Fig. 7. At low shock velocities, 
the composite profile is virtually coincident with the Maxwellian of the protons and so, to a 
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very good approximation, Av{vs, 0) = 3i's/4. However, for Vg > 2000 km s ^, the composite 
profiles (Fig. 2) begin to shift noticeably toward the narrow line, so that A.v{vs, 0) < 3vs/4:. 

It should be possible to infer the aspect angle of the shock by comparing the mea- 
sured FWHM with A^;. For example, KWC87 measured FWHM = 1800 ± 100 km s"^ and 
Av{vs, 9y) — 238 ± 18 km s~^ in a filament in Tycho's remnant. From this result and Fig. 7, 

wc can infer Vs = 2409 ± 156 km s"^ and 9^ = SA^s^^^o for Model F, consistent with KWC87's 
estimation of 90° — 6^ ^ 6°. Of course, the observations will tend to select bright filaments 
that are viewed nearly edge-on. 

The measured ratio of broad to narrow line Ha intensities, 3?6,„(Hq;) = Ib(lia) / Ifj(lia) , 
is listed in column 5 of Table 1. As first noted by CKR80, this ratio is sensitive to the optical 
depth in the in the Lyman series. In choosing the cascade matrices in equations (20) and 
(22), we consider two limits. In Case A, all Lyman transitions are optically thin. In Case B, 
they are optically thick; therefore, any Ly/3 photon that is emitted will be trapped until it 
splits into LycK-t-Ho;. Likewise, Ly7 must be converted to Lja+B.p or Lya+Ha+Pa, and so 
on. In Balmer-dominated shocks, the narrow line emissivity of Ha for Case B is increased 
by a factor 2.5 compared to Case A, and that of H/3 by a factor ^ 1.9. GOl modeled the 
ionization structure of the shock fronts in three SNRs and performed Monte Carlo simulations 
of Lyf3 and Ly7 trapping. They find that most of the broad Lyman photons escape from the 
shock, but that narrow Lyman photons are likely to be absorbed and converted into Balmer 
lines both behind and ahead of the shock. Therefore, for most SNRs (except SNR 1987A), 
we are making a good approximation to adopt Case A for the broad lines and Case B for the 
narrow lines. For a theoretical upper limit on we assume Case A for both components. 

In Figure 8, we plot the measurements of 3?6,„(Ha) and Vs inferred from the measured 
FWHM for several SNRs. We do this for both the non-equilibrated Model N and the fully 
equilibrated Model F. As in columns 6 and 7 of Tabic 1, the uncertainties in the inferred Vg 
are due entirely to the measurement errors. We also plot, for both models, the relationship 
3?5,„(Hq;) vs. Vs predicted by our theory, as well as the same relationship from equation (4) 
of S91 (using our rate coefficients). 



which follows from the model originally presented by CKR80. The Case A and B efficiencies 
{ea and es) are given in CR78, while we adopt the mean value of the fraction of charge 
transfers to excited states yielding Ha photons to be ~ 0.03. In both cases, we assume 
Case A for and Case B for /„. 

One can see from Figure 8 that the data fit the models fairly well, and that our derived 
relationship of 3?;,,n(Ha) vs. Vg does not differ very much from that derived by CKR. In the 
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cases of SN 1006, Kepler and RCW 86 and 0548 — 70.4, the fits do not discriminate very well 
between Models N and F, though Model F provides a better fit to the Cygnus Loop than 
Model N. The remnants of Tycho and 0519-69.0 are marginally fitted by Models F and N, 
respectively. With the exception of SN 1006, we also note that the data cannot be fit with a 
model in which Case A applies to /„ as well as lb- In that case, the theoretical curves would 
be elevated by a factor pa 2.5 (Fig. 8), and would not fit the data at all. 

We can make a similar comparison of theory to data with the H/9 line. The corresponding 
fraction of charge transfer to excited states yielding H/5 is estimated to be Q/s ^ 0.01 (Bates 
& Dalgarno 1953). Figure 9 shows both the measured and theoretical values of 3?f,,ra(H/?) vs. 
Vs- Again, the difference between our model and that of CKR80 is not great. The fit of the 
theory to the data does not discriminate between Models N and F, except for the Cygnus 
Loop where the latter is favored. 

Note that the ratios ^b,n{iic() and Kb,„(H/3) are fairly insensitive to shock velocity in 
the range 200 < t's < 2000 km s~^. However, for Vg > 2000 km s~^, these ratios do become 
quite sensitive to Vg. It is worthwhile to note that the Monte Carlo simulations of GOl and 
G02, which model Lyman line trapping, result in tighter constraints on ^b,n- 

5.2. SNR 1987A 

SNR 1987A has a double shock structure, consisting of a forward and a reverse shock. 
In contrast to other SNRs, the Balmer-dominated shock in SNR 1987A comes, not from 
nearly stationary neutral hydrogen overtaken by the supernova blast wave, but from freely- 
expanding neutral hydrogen atoms in the supernova debris crossing the reverse shock. The 
debris crosses the reverse shock at velocities ~ 12,000 km s"-*^ (M03; H06). In the frame 
of the observer, the forward shock or blast wave is moving at between 3500 and 5200 km 
s'^ (Manchester ct al. 2002; Park ct al. 2002; Michael et al. 2002). The reverse shock is 
moving at about 80% of the velocity of the forward shock (M03). Hence, in the frame of 
the observer, the post-shock ions are relatively slow particles moving at velocities of between 
5100 and 6120 km s^^. The "shock velocity" in this case is interpreted as the velocity of the 
freely-streaming hydrogen atoms in the rest frame of the reverse shock, i.e., 7840 < Vg < 9200 
km s~^. 

In the terminology of this paper, the "narrow" line emission resulting from excitation 
of hydrogen atoms initially crossing the reverse shock will have the greatest Doppler shift, 
corresponding to the free-streaming debris. There is a unique mapping of the Doppler shift 
of this narrow emission to the depth of the reverse shock measured from the midplane of the 
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supernova debris along the line of sight. Using the Space Telescope Imaging Spectrograph 
(STIS), M03 and HOG (and references therein) mapped streaks of high velocity emission 
in both Ha and Lya in SNR 1987A, which H06 called "surface emission". This surface 
emission is evidently equivalent to narrow line emission from atoms crossing the reverse shock. 
H06 also detected Ha and Lya emission from the same location having substantially lower 
Doppler shifts than the surface emission. If that emission came from freely-expanding atoms 
in the supernova debris, those atoms would have to reside beneath the reverse shock surface. 
Accordingly, H06 called this emission "interior emission" . That was an unfortunate choice 
of terminology. At least part if not all of this emission must come from more slowly moving 
hydrogen atoms immediately outside the reverse shock surface, which result from charge 
exchange reactions of hydrogen atoms in the supernova debris with protons in the shocked 
gas. This "interior emission" is analogous to the broad line emission in the terminology of 
this paper. 

Following M03, we have adopted a model for SNR 1987A having xhc = nHe/nn = 0.2. 
As before, we consider two variants. Model N {Te/Tp = 0.25, feq ~ 0.33) and the fully 
equilibrated Model F {T^/Tp = 1). In SNR 1987A, the debris crossing the reverse shock has 
a Sobolev optical depth Tiya ~ 1000 (M03). Therefore, we assume that Case B is a good 
approximation, not only for the narrow line emission, but also for the broad line emission. 

Since the broad and the narrow line emission from SNR1987A comes from the expanding 
debris as it crosses the reverse shock, the line profiles will depend on the hydrodynamics of 
the reverse shock surface, which we do not consider here. However, one can measure the 
ratios of broad ("interior") to narrow ("surface") emission and of Ha to Lya emission (H06). 
By comparing those ratios to those predicted by our theory, we can test the hypothesis that 
the interior emission is indeed the result of charge transfer reactions between hydrogen atoms 
and protons at the reverse shock. 

Figure 10 shows these ratios as functions of shock velocity for our model of SNR 1987A. 
The sohd curves show the photon emission ratio ^hl — /(Ha)//(Lya). We see that '^hl ~ 
0.19, a result that is almost independent of Vg and almost the same for interior and surface 
emission, both for Model N and Model F. This value is consistent with the value ^hl ~ 0.2 
estimated by CKR80 and ?R:hl ~ 0.21 estimated by M03. 

The dotted curves in Fig. 10 show the number, Niya, of Lya photons (both surface 
and interior) emitted per hydrogen atom crossing the reverse shock. In the rest frame of the 
reverse shock, these atoms cross the shock with a velocity 7840 < f ^ < 9200 km s~^. For Vg 
in this range, N^ya ~ 1.1 for both Models N and F. These values may be compared to the 
value NLya ~ 1 estimated by M03. 
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The dashed curves in Fig. 10 show the ratios 3?b,n(HQ;) of interior to surface emission. 
For the range of shock velocity expected in SNR 1987A, the curves predict ^b^nO^ct) ~ 0.08 
for both Models N and F. These values are more than an order of magnitude smaller than 
the ratio K(,^„(Hq!) ~ 1 observed by H06. 

6. DISCUSSION 

For most SNRs, the models presented here agree fairly well with previous models and 
seem to account for the observed data of both broad hue FWHM and ratio 3f?fe,n(HQ;) fsi 1 of 
broad-to-narrow line intensities. However, these same models fail to account for the observed 
ratio 3?6 „(Ha) of interior-to-surface emission for the conditions of the reverse shock in SNR 
1987A. Either the model is wrong, or the interior Ha emission from SNR 1987A is dominated 
by some excitation mechanism other than charge transfer at the reverse shock. 

As we mentioned in §2, the most questionable assumption of the model is that the 
hydrogen atoms enter a proton plasma which has been suddenly decelerated according to 
shock jump conditions. This assumption is flawed — the approximation is good only when 
the neutral fraction is small (see Lim & Raga 1996 for a discussion on the effects of larger 
neutral fractions). Indeed, the reverse shock has a transition zone of finite thickness Az ~ 
UpUj^p, within which the ionization, excitation and charge transfer reactions occur. Consider 
the first hydrogen atoms to enter this zone from upstream. When they become ionized, the 
resulting protons must have streaming velocities comparable to that of the hydrogen atoms. 
Within the transition zone, the newly created protons must be decelerated and compressed 
so that they emerge on the downstream side of the zone with velocity Vp = Vs/4:. But most 
of the ionizations must occur when the relative streaming velocity of the hydrogen atoms 
and the protons is substantially less than the relative velocity of 3^3/4 that is assumed in 
the present model. Given the fact that 3?6,n(HQ;) decreases precipitously for Vg > 2000 km 
s~^ (Fig. 10), we see that our model might grossly underestimate the actual value of this 
ratio in the case of SNR 1987A. 

To properly interpret the Ha and Lja emission from SNR 1987A, we must construct a 
model for this transition zone similar to that of Whitney & Skalafuris (1963) and Skalafuris 
(1965, 1968 and 1969). In a sense, the analysis presented here is a warm-up to that task. 
Once we have such a model, we can use the results of the present analysis to calculate the 
emissivity and profiles of broad and narrow lines throughout the transition zone. 

The assumption that a Balmer-dominated shock can be treated as a discontinuous jump 
is faulty for all SNRs, not just SNR 1987A. Accordingly, we may ask why the model seems 
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to work for most SNRs, for which Vg < 2000 km s~^. If one allowed for the transition zone in 
constructing a model for actual SNRs, one would construct a weighted average of $R5 „(Hq;) 
over relative velocities ranging up to 3vs/4:. But, as Fig. 8 shows, the ratio 3?f,^„(HQ;) does 
not vary greatly for 200 ^ < 2000 km s~^. Therefore, the value of 3?;,^„(Hq;) appropriately 
averaged over a range of relative velocities would still agree fairly well with the observed 
values. 

On the other hand, these considerations bring into question the validity of the relation- 
ship between the broad line FWHM and the shock velocity predicted by the current models. 
In a more realistic model that takes into account the finite thickness of the transition zone, a 
greater fraction of the Ha emission would take place where the relative velocity is less than 
3vs/4:. The result would be a reduced value of the FWHM of the broad emission for a given 
value of Vs — or, equivalently, a greater value of Vg required to account for a given FWHM. 
This consideration suggests that the current models may be under-estimating significantly 
the shock velocities of SNRs. 

The present model is also inadequate to interpret current observations of the ratio 
^b^niLya) in SNR 1987A. The observed brightness of Lyct is modified by resonant scattering 
in the supernova debris. The resonance scattering can double the brightness of Lya on the 
near side of the supernova debris and greatly suppress it on the far side, as pointed out by 
HOG. However, since the Lya photons are produced in the transition zone, and their mean 
free paths are less than Az, we need a model of the transition zone to account properly for 
this resonant scattering. 
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expertise of James McKown and Peter Ruprecht, who provided impeccable technical support at JILA, without which this piece 
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Consider a hydrogen atom in the distribution function /q. If it undergoes an infinite 
number of charge transfers, the average number of excitations it will have is 



A. 



APPENDIX 



A.l. ADDITIONAL EQUATIONS FOR THE REACTION TREE 
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where the second term accounts for excitations after one or more charge transfers. In the 
case of a reaction tree with no 7] reactions, equation (Al) simply reduces to the famihar 
Peo/{1 — Peo)- Similarly, the average number of ionizations and charge transfers per atom 
are 



and 



OO J 



T 



j=0 k=0 



respectively. 

The total rate is obtained in an analogous way: 



R 



■X,total 



If we make the approximation that Rx^ — Rx for i > 1, then the preceding expression 
reduces to: 

1 Pto 



Rx, total = i_ [P^o + Pxj ■ (A5) 

The zeroth order rate is then Rxq/ (1 — -Peq)- Examples of Rxa/nni)- — Peq) and Rx,totai/nH 
are plotted in Fig. 11. 



A.2. CASCADE MATRIX 

Consider the coUisional excitation of a hydrogen atom from the ground state to the 
level n. Let P„„' be the probability that the population of n is followed by a direct radiative 
transition to n' (Seaton 1959), given by 

Pun' = (A6) 

where 

n-l 
i=no 

and we have Uq—I for Case A and no—2 for Case B. 
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Let Cnn' be the probability that the population of n is followed by a transition to n' via 
all possible cascade routes. Following Seaton (1959), the cascade matrix elements are 

n-l 

Cnn' = PniCjn', (-^8) 

i=n' 

where, for any a >b, Cab=^ if a = 6. An analytically more illuminating form is 

rt— 1 n— 1 n— 2 

Cnn' Pnn' ~l~ ^ ^ Pni]_Piin' ~l~ ^ ^ ^ ^ Pnii Piii2 Pi2n 
ii=n'+l ii=n' +2 i2=n' +1 

n—1 n—3 n+3 rt+2 

"I" ^ ] ^ ] ••• ^ ] ^ ] Pni\Pni2---PiN-iiNPiNn' 

ji=n— 2i2=n— 2 i/if_i=n+2 jjv=n+l 

+-Pn,n-l-Pn-l,ra-2----Pn+2,n+l-Pn+l,n) (^9) 

where for any a> b, Pa,b = -Fob and Pab—^ for a = b. 

Determining the cascade matrix elements requires knowing the Einstein A-coefficients. 
In the electric dipole approximation, they are (Pengelly 1964) 

^ n— 1 

= + l)Anl,n'l', (AlO) 

where the quantity Ani,n'i' is given by (e.g., Pengelly 1964; Brocklehurst 1971) 

A„,,„.^ 2.6774 X 10»(± 4)' !=|fi^|<f|^. (All) 

An extensive literature exists on how to calculate the radial integral, R^j . Gordon (1929) 
first provided exact expressions for it, in terms of hypergeometric functions. Multiple authors 
have since developed different ways of computing this formula, mostly in tabular form (e.g.. 
Green, Rush & Chandler 1957; Goldwirc 1968; Menzel 1969; Khandelwal & Fitchard 1972). 
Malik, Malik & Varma (1991) rewrote Gordon's formulae in terms of associated Laguerre 
polynomials, yielding 

^n' ; 1 ^ ^ ^J-^ f n + l\ f n' -\-l — 1 \ + Z — A — 3\ ^ ^ ^ ^, 

K/-' = ^-.^-.E ( A j(„'-„ + A + 2)( 21-2 j»^m + r, + r3 + r, + r,) 

A=Ao 

(A12) 

and 

Rn,i =^^^3 E A jU-- + A + 2jl 21 jMti+t, + t3 + t, + T,). 

A=Ao 

(A13) 
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where Aq = max (0, n — n' — 2) and (^) = a\/{a — b)\ b\ for any a, 6 e Z is the binomial 
coefficient. 



Tfie various functions used in equations (A12) and (A13) are 
Fi = 



and 



n + I 
21-2 



n' + 1-1 
21-2 



l + l){n- 


l){n' 






nn' 


n+l ^ 


_{n + n'y _ 





-1/2 



n — n 
n' + n 



n' + l + f n + l 



21 



{n'-l + l)(n' - l)(n' - I - l)(n - 1) 



-1/2 



21 

Ti = -y\2 + + l)(j + / - 2)(j - / - 1), 
n = 2y{j + l){j + l- l)(2j - / - l)(n' - J + 2), 
Ts = -6j(i + l){n' - j + l){n' -j + 2), 
T, = -{2j + l + l){n'- j){n' - j + l){n' -j + 2), 

y 

n = -\{n' - J - - j){n' - j + l)(n' -j + 2), 

y2 



ti = -y'U + -l-m-l- 2)0' -I -I). 

t2 = 2y{n' -j + 2){j-l- - I - 2){2j + I), 
ts = -6i(n' - j + l)(n' - i + 2)(i - / - 1), 

U = -(n' - j){n' - J + l)(n' - j + 2)(2j - /), 

y 

(n — n'Y 

Ann' ' 
n' — n 

~ 2n ' 
j ^n-\. 



y 



(A14) 



(A15) 



where wc note that the last term in Fi was erroneously written as / in Malik, Malik & 
Varma (1991). With the correct version, equation (A12) reduces to the following, exact 
form (Hoang-Binh 1990): 



Rl-'fT' = V(2n-l)(2n-2) 



n{n — 1) 
2n- 1 



4n(n — 1) 
(2n - 1)2 



(A16) 
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Finally, we display the resulting cascade matrices in (Fig. 12). We have assumed that 
the excitations are statistically distributed among the angular momentum states. We note 
that this is a reasonable assumption for charge transfer into excited states, but not necessarily 
so for impact excitation, which favors np states, especially at high energies. 
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Table 1: FWHM of Broad Component of Hct, h/In and Predicted Shock Velocities 



SNR 


Reference 


FWHM 


Vs 


cited^ 


h/In cited 


Vs (Model N)^ 


Vs (Model F)^ 






(km s~^) 


(km s ^) 




(km s~^) 


(km s~^) 


SN 1006 


G02 


2290 lb 80 


2865 


— 3580 


0.84i°:°? 


2509 ± 111 


2981 ± 133 


SN 1006 


S91 


2310 ± 210 


2400 


— 3240 


0.73 ± 0.06 


2537 ± 292 


3014 ± 354 


SN 1006 


KWC87 


2600 ± 100 


2800 


— 3870 


0.77 ± 0.08 


2940 ± 149 


3545 ± 183 


Kepler 


BLV91 


1500 


1550 


— 2000 


0.72 ± 0.37 


1528 


1809 


Kepler 


F89 


1750 ± 200 


1670 


— 2800 


1.1 ± 0.25 


1806 ± 226 


2154 ± 294 


Tycho 


GOl 


1765 ± 110 


1940 


— 3010 


0.67 ± 0.1 


1823 ± 122 


2177 ± 169 


Tycho 


S91 


1900 ± 300 


1850 


— 2500 


0.77 ± 0.09 


1973 ± 372 


2383 ± 453 


Tycho 


KWC87 


1800 ± 100 


1930 


— 2670 


1.08 ± 0.16 


1862 ± 111 


2230 ± 153 


RCW 86 


GOl 


562 ± 18 


545 


— 793 


1.18 ± 0.03 


536 ± 18 


641 ± 21 


Cygnus 


GOl 


262 ± 32 


235 


— 395 


0.59 ± 0.3 


247 ± 30 


296 ± 36 


0505—67.9 


S91 


580 ± 70 


480 


— 640 


>0.7 


554 ± 69 


662 ± 82 


0519—69.0 


S91 


1300 ± 200 


1180 


— 1580 


0.8 ± 0.2 


1305 ± 223 


1554 ± 255 


0518 70.1 


S91 


760 ± 1 10 


670 


890 


1.1 ± 0.2 


732 ± 1 10 


871 ± 168 



Range of shock velocities quoted for zero and full electron-ion equilibration. 
Derived for edge-wise shock fronts. 



100 1000 10000 

Relative velocity (km s"') 



Fig. 1. — Cross sections for impact excitation (E), ionization (7) and charge transfer (T) for 
interactions between neutral hydrogen atoms and various species. The subscripts "p", "e" 
and "a" refer to proton-atom, electron-atom and alpha particle-atom reactions, respectively. 
The E, I and T (n = 1 only) data are from Janev & Smith (1993), while that for charge 
transfer to the n = 2 state (T*) are from Barnett (1990). 
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Fig. 2. — Face-on projected profiles, 0, for the proton and hydrogen atom velocity distri- 
butions, using arbitrary relative units. The protons are in Maxwellian distributions (0Af), 
while (f)i and (f)comp are the atomic distributions after i and one or more charge transfers (i.e., 
composite), respectively. The velocity normal to the plane of the shock front is denoted by 

Vz- 
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Fig. 3. — Reaction rate coefficients for impact excitation (Ei), ionization (Jj) and charge 
transfer (Tj), for various species and both Models N (top row) and F (bottom row). The 
subscripts "p", "a" and "e" refer to proton-atom, alpha particle-atom and electron-atom 
reactions, respectively. Shown are the rate coefficients for reactions involving atoms from 
the original beam (i = 0) and those which are the result of one charge transfer reaction 
{i = 1). 
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^EjJiJi (i>i) 

►EiJiJi (i>2) 

I 

■►EiJiJi (i>2) 
■►EiJiJi (i>3) 



Fig. 4. — Reaction tree of impact excitation (Ei), ionization (Jj) and charge transfer (Tj). 
The index i denotes the number of times a hydrogen atom has experienced charge transfer 
reactions. 




Fig. 5. — Reaction probabilities, P^^, Pe^ and P/., for excitation, ionization, and charge 
transfer, respectively. Shown are the probabilities of the reactions involving the hydrogen 
atom beam and the mean probabilities (see text). 




Fig. 6. — Average number of excitations (Ne), ionizations (Nj) and charge transfers (Nt) 
per atom. The zeroth order quantities are Neq — Peo/(1 — Peq), -^/o — Pio/i^ ~ -^Eo) 
Nt, = PtJ{1-Pe,). 
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Fig. 7. — Full width at half-maximum (FWHM) of the projected profiles of fcomp- These 
are shown for both face-on and edge-wise shock fronts. Also shown are the velocity offsets 
between the peaks of the narrow and the face-on broad components {(pcomp)- The data 
points represent velocity offsets measured for Tycho by Kirshncr, Winkler & Chevalier (1987; 
238 ±18 km s^"*^) and Ghavamian ct al. (2001; 132 ±35 km s^^), which are denoted by {Kj^f, 
Kp) and (Gat, Gp)-, respectively; the subscripts stand for Models N and F. We assume a 
face-on shock front to derive the shock velocity range for each data point. 
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Fig. 8. — Ratio of the broad to narrow Ha emission, ^^^^(Ha), vs. shock velocity, Vg- Each 
pair of thick curves represents our lower and upper bounds for K^^^ (denoted by the label 
"HM06"). The "CKR80" curve is obtained using the expression of Chevalier, Kirshner & 
Raymond (1980) and Smith et al. (1991). The data represent SN 1006 (Ghavamian et al. 
2002), Kepler (Fesen et al. 1989), 0519-69.0, 0548-70.4 (Smith et al. 1991), lycho, RCW 86 
and Cygnus (Ghavamian et al. 2001), assuming edge- wise shocks. 




Fig. 9. — Same as Fig. 8, but for B.(3 emission. The data represent SN 1006 (Ghavamian et 
al. 2002), Tycho, RCW 86 and Cygnus (Ghavamian et aL 2001), again assuming edge- wise 
shocks. 
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Fig. 10. — Various ratios of the Lya to Ha and Lya to ionization rate coefficients for 
conditions relevant to SNR 1987A (Case B and 20% alpha particles by number compared to 
hydrogen). "Total" refers to the sum of surface and interior emission, as defined by Heng at 
al. (2006; see text). 




Fig. 11. — The total rate coefficient, Rx,totailtT-Hi weiglited by tfie fractional abundance of 
each species (protons, electrons and alpha particles), for excitation, ionization and charge 
transfer. The zeroth order counterparts (see §A.l) are denoted by "0th". 



-35- 



1.00 F 



0.10 



0.01 



1.00 



0.10 



0.01 



n' = 2 
n' = 3 
n' = 4 
n' = 5 



n' = 2 
n' = 3 
n'=4 
n' = 5 



15 



20 



Fig. 12. — Cascade matrices for Case A (top) and B (bottom). Each matrix element, Cnn', is 
the probabihty that the population of n is followed by a transition to n' < n via all possible 
cascade routes. 



